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ABSTRACT 


Ths effects of environmental perturbations on the attitude of a 
slow tumbling earch-oriented satellite are investigated. The environ- 
mental perturbations considered were aerodynamic drag, gravity-gradient, 
solar radiation pressure, and magnetic torques. The Euler attitude 
equations were solved numerically for the Skylab spacecraft. Results 
are presented for both torque-free motion and for cases in which 
aerodynamic and gravity-gradient torques are acting in a slow tumble 
mode. Simulations show gravity-gradient effects on satellite momentum 
to be cyclic and to increase the precession rate of the angular momentum 
vector about the radius vector. This also tends to align the minor axis 
along the radius vector. Aerodynamic drag initially decreases angular 
momentum, slowly precesses the momentum vector about the radius vector, 
and finally drives the satellite into an unstable mode. Combined 
gravity-gradient and aerodynamic torques reduce angular momentum and 
energy, and induce a steady precession rate of the mpmentum vector about 


the radius vector. 


CHAPTER I 


INTRODUCTION 

With the coining of the space shuttle the opportunity for retrieving 
and repairing satellites will become feasible. With retrieval capability, 
future satellites may be designed with docking ports to be used by a 
retrieval device for attaching to the satellite. To design these 
devices and the location of docking ports on the satellite, the final 
attitude state of the passive satellite must be determined. 

The research in this thesis involves a numerical study of the 
general attitude motion for an asymmetric satellite due to environmental 
perturbations. These perturbations included gravity-gradient, aerodynamic 
drag, solar radiation pressure, and magnetic torques. 

The general analytic solution to Euler's moment equations are 
elliptical functions for the torque free case.^^^ General solutions to 
Euler's equations with complicated torque functions do not exist and 
solutions are primarily numerical. Some solutions have been obtained by 
linearization for special cases. Here the equations could not be 
linearized to study detailed motion. However, Euler's moment equations 
were solved numerically using a fourth order Runge-Kutta and an Adams- 
Bashforth predictor-corrector integration technique. The dynamical 
state of interest was a slow tumble mode, defined as small angular rates 
about the three body axes. Initial conditions are explained in Appendix 
A. The satellite was considered to be a rigid body with no control 
system functioning. 

The asymmetric satellite studied in this research was the Skylab 
spacecraft. Its attitude control system is assumed to be shut down. 


Linearization of the equation of motion was possible in order to study 

C 2 3 A ^ 

stability under effects of gravity-gradient and aerodynamic drag* * ’ ^ 
This stability analysis has indicated the spacecraft to be unstable in 
the presence of gravity-gradient with aerodynamic torques. 


CHAPTER II 


COORDINATE SYSTE^^S 


This section describes the coordinate systems used in determining 
the position and attitude of the satellite. The inertial coordinate 
system is shovn in Figure 1 with its origin at the Earth’s center. 
lies in the equatorial plane pointing in a positive direction away from 
the center of the earth along the vernal equinox. is perpendicular 
to the equatorial plane and positive northward. Xj. is in a direction 
forming a positive right handed coordinate system. The inertial 
coordinate system will be denoted as 



( 1 ) 


The relationship of the apparent motion of the sun about the earth 
is shown in Figure 2. The unit vector from the earth to the sun in the 
geocentric inertial coordinate system is 


L = sin 9 cos i i_ + sin 9 sin i j, + cos 9 k_ (2) 

s s si s s i si 

where 


- -360 n 
*^3 365.24 AE 


(3) 


is the number of days after vernal equinox, i^ is the inclination 
of the ecliptic plane to the equatorial plane of the earth (obliquity 
of the ecliptic). 

The relationship between the orbital and inertial coordinate 


system is shown in Figure 3. The transformation from inertial to 
orbital coordinate system is defined by three rotations. 







7 


Y' 



Figure 3. Inertial and Orbital Coordinate System (continued). 
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(4-B) 
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The orbital coordinate system will be denoted by: 


'"'lo' 


X 


(5) 


From Figure 3, the following transformation exists: 
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Z’ 
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( 6 ) 


From transformations (A) and (6) the transformation from inertial to 
the orbital coordinate system becomes 


[L] 


o->I 


(cosocosicosg-sinasinB) (cososini) (-cosacosising-sinacosg) 

(sinicosg) (-cosi) (-sinisinB) 

(-sinacosicosS-cososinB) (-sinosini) (sinacosisinB-cosacosB) 

(7) 


where 
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B is the right ascension of the orbit ascending node. The nodal 
regression rate is assumed constant, 
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B * - 0.001637 (^)‘ 
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1800P 


cos i 


(9) 


10 


where P is the satellite orbital period in hours, i is the satellite 
orbit inclination with respect to the equatorial plane, and a is the 
satellite true anomaly; its initial displacement is measured from the 
ascending node for a circular orbit. 

The body coordinate system is located at the center of mass of 
the satellite along its principal moments of inertia. The rotation 
from the orbital to the body coordinate system is illustrated in 
Figure 4. The order of transformation is defined as 


X 
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The three Euler rotations are 


A. Rotation about the Z axis 
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cos ^ sin il) 

-sin {jj cos ijj 

0 0 



( 10 ) 



(10-A) 
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B. Rotation about the axis 
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C. Rotation about the X 2 axis 
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The transformation from orbital to the body coordinate system can be 
expressed as 


[L] 


b-+o 


(cos0cosi|>) (cos6sini|j) (-sin0) 

(sincf)sin0cosi|J-cos(|)sinij^) (sin(|’sin0sin^cos({>cosi[') (sin({)cos0) 
(cos(|)sin0cosijH-sin(J>sim|;) (cos4>sin0sini|;-sin4)cosi|O (cos(})cos0) 


( 11 ) 


where 




— — 


, 

X 

0 



Y 

0 



Z 

0 


(11-A) 


The above transformations are basic in determining the attitude 
and position of the satellite. All rotations are defined positive by 
the right handed rotation rule. Any other transformations required are 
developed as needed. 


REPRODUCIBILITY OF THE 
ORIGD'JAL PAGE IS P(X)R 
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CHAPTER III 
DYNAMIC ANALYSIS 


The Liotion of a satellite about its center of mass is described 
by Euler's moment equations which for the principal axes are^^^ 

Aw + u) (0 (C-B) = M 

X y z X 

B (I) + 0) u) (A-C) = M (12) 

y X z y 

C 0) + 0) W (B-A) = M 

z X y z 

A, B, and C are the principal moments of inertia about the x, y, and z 

body axes, respectively, w , to , and to are the satellite angular 

X y z 

velocities. M^, and are the perturbing moments such as gravity- 
gradient torque. 

• • • 

Euler angular rates tjJ, 0, cj) describe the motion of the satellite 

with respect to a reference coordinate system i^hich in this report is 

the orbital coordinate system. Euler angular rates can be expressed as 

a function of the body angular velocities and the Euler angles. 

d) = 'CO + (co sind) + (O coscb) tan 9 
X y z 

6 “ “ (O^sincj) (13) 

= (cOySincj) + (O^costt)) sec 9 

The body angular velocities can also be written in terms of Euler rates. 
(0 = 4> - di sin 0 


cOy = 0 coscj) + d) cos 0 sin 4> 


(0^ “ COS0 cosc}) - 0 sin (j) 
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Equation (12) can be expressed as 
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(15) 


Knowing the perturbing moments as a function of the Euler angles 
equations (12) and (13) can be solved numerically. These equations 
describe the attitude motion and orientation of an asymmetric satellite. 

The angular momentum, h and energy, T of the satellite are 
computed from 

h = Ao3 i, + Bo3 j- + Co3 k, (16) 

X b y b ZD 



(17) 


The position of the angular momentum vector in the orbital coordinate 
system with respect to the orbit radius vector and velocity vector 
tangent is shown in Figure 5. The angle 6 is the angle between the 
radius vector of the earth and the satellite momentum vector. 

1 h 

6 = cos"-*- (^) 0° < 5 < 184>° (18-A) 

n O 

h^ is the Z angular momentum component in the orbital coordinate system. 
The angle X defines the angle between the angular momentum projection 
onto the x^, y^ plane and the orbital velocity vector tangent. 


Orbit Normal 



Figure 5. Nutation and Precession Angles in the 


Orbit Coordinate System. 
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X = tan"^ (^) 0° < X 360° (18-B) 

y o 

and are the x and y angular momentum components in the orbital 
coordinate system, respectively. 5 and X will be noted as the nutation 
and precession angles in the orbital coordinate system. 

The position of the minor axis with respect to the radius vector 
will be noted as the nadir angle, n* The nadir angle is defined as 

n = f + 0 


(19) 
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CHAPTER IV 

ENVIRONMENTAL PERTURBATION MODELS 


Gravity-Gradient 

Gravity-gradient torque is one of the major environmental 
perturbations for asymmetric satellites in near earth orbits. The 
gravity-gradient effect is a function of altitude, mass distribution, 
and the satellite orientation. Here the gravity torque model assumes a 
spherical earth neglecting anomalies due to its asymmetric mass 
distribution. 

Under the assumption of a spherical earth, the gravity force 
field can be expressed as 



r 

and the gravity-gradient torque as 


( 20 ) 


(dM)g » p X dPg (21) 

Figure 6 illustrates the coordinate system used in this derivation, p 
is the vector from the satellite center of mass to the mass element dm. 

p = xi^ + yj^ -f (22) 

Subscript b indicates the satellite body coordinate system, r is the 
satellite radius vector from the geocenter to the mass element dm. 

From Figure 6 the satellite radius vector R can be written as 
R = -R 

o 

Transforming equation (23) to the satellite body coordinate system 
using equation (11) , the radius vector can be expressed as 


( 23 ) 


\Mass Element (dm) 



Figure 6. Coordinate System Used in Gravity-Gradient 

Torque Derivations. 
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, r-. ' -1.'?. : 




R^ = R sin6 - R sinc|) cos9 I|j ” ^ cos(f> cos9 


Noting that r =* R + p equation (21) becomes 


(dM)^ = p X [ - -Hj (R + p) ] 
r 


With the following approximations; 


== r2 [ 1 + 2R^ ] 

R 


1- ~ 1. r 1 + 2R-B T -3/2 


1 r 1 3R»p , 

3 “ 3 ^ 2 ^ 

r"* R-^ R^ 


p/R « 1 

Equation (25) becomes 


dM, = ^ [ 1 - [p X (R + p) ] dm 


° r2 


Since the coordinate system is located at the satellite center of mass 
the products of inertia are zero 


p dm S' 0 


I = I = I = 0 
xy xz yz 
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Integrating equation (30) results in the following equation for the 
gravity-gradient torque 


M, 


31L. 

2R^ 


[sin2(|) cos*^6 (C-B)i^ + sin29 cos<J> (C-A)j^ 


+ sin20 sin^) (A-B)^ ] 


(33) 


Equation (33) describes the gravity-gradient effect in satellite body 
coordinate system. 

The gravity-gradient torque in the orbital coordinate system can 
be found by transfoirming equation (33) to the orbital coordinate system 
using the inverse of equation (11). The resulting torque equations in 
the orbital coordinate system are 


(ra ) =* [(A-B)sin29 slni|J + (B-C) {sin29 sinij; cos^(J) + sin2<}) cos9 cosi^}] 

X o 2R-^ 

(m ) =[(B-A) sin29 cosii; +(C-B) {cos 6 sin29 cosi|> - ,sin2(j) cos9 sinij^}] — r 

y ° 2R^ 


(m ) 

z o 


0 



(34) 


(iti ) j (m ) , and (m ) are the x, 
'• X o y o z o 

in the orbital coordinate system. 


y, and z gravity-gradient components 


Aerodynamic Drag 

Aerodynamic drag is a major perturbation for near earth satellites 
with altitudes of 800 km or less. Drag is a function of atmospheric 
density, angle of attack, satellite velocity and satellite shape. For 
complex satellite structures the satellite is divided into components of 



= i S Pa '^V 

is atmospheric density. is the aerodynamic drag coefficient. V 
is the velocity at the surface element relative to incident stream. Yy 
is the angle of attack of element dS. Equation (36) is the basic 
aerodynamic force equation. Reference (8) discusses aerodynamic force 
using normal and tangential momentum transfer coefficients to replace the 
drag coefficient. 

(9) 

For this report the aerodynamic drag model was developed by NASA 

for the Sky lab vehicle. The drag model for Sky lab was derived from data 

based on free molecular flow theory with a Knudsen number greater than 

10. Three drag moment coefficients (c , c , c ) as depicted in Figure 7 

y z 
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were computed for a number of orientations as a function of the angles 

(a . d> ) . A Fourier series curve fit formula for c , c , c was 
a a X y z 

derived from the above data as a function of the angles ot and (j> , 

'3- H 

c . c . and c are the roll moment coefficient, pitching moment 
X* y z 

coefficient, and yaw moment coefficient, respectively. Ci and <j) are 
the angle of attack and roll angle, respectively, as defined in Figure 7. 
The resulting Fourier drag coefficient equations are 


A ((j) ) m 

c(a B.(tj)„)sin ia ] 

a a 2 1=1 1 a a i a a 


(37) 


where 


““^♦a '=alj 


(38-A) 


(for j = 2, 4, 6, and i = 1, 3) 




\lo 

2 


n 

+ Z 

j=l 


(for j 


[a... cosj(}) +b sinjcf) ] 
bij a bij a 

= 1, 3, 5, and i = 1, 3) 


(38-B) 


a and b are coefficients from Appendix B. The vehicle moment 

Ij ij 

equations are computed from 


M = c q A - D - 
X X ref ref 


M = cqA -D p 
V y ref ref 


M = c q A . D . 
z z ref ref 


where 


q = I pv^ 


( 39 ) 


A ^ and D ^ are reference Area and Diameter of Sky lab respectively, 
ref ref 

and are listed in Appendix C. 

The roll angle (ji and angle of attack a were computed from 

Si 3 

^ V 

d) = tan~-^ 0 < 4) < 360° (40) 

a V fa a 

z 

and 

- V 

a = cos"-^ 0 < a ^ 180° (41) 

a V D a 

as shown in Figure 7. Assuming a circular orbit for Sky lab, where V 

o 

is the orbital velocity, 


V = V i 
o o o 

Expressing the velocity in body coordinates results in V^, V^, and 

V as required in equations (40) and (41) 
z 


V 

X 

V 

y 


V 

z 



cos6 cosif) i^ 
[sin4> sin6 cos^^ 
[cos(j) sin0 cosi|; 


- cos4> sinilJ] j, 


+ sintj) sinijj] k^ 


(42) 


Appendix B lists the aerodynamic coefficients from Reference (9) 
used in computing the aerodynamic drag torques. The drag coefficients 
in Appendix B are based on Skylab with the auxiliary thermal shield, ATM 
solar arrays, and orbital workshop solar panel No. 1 deployed as shown 
in Figure 8, 

The atmospheric density was calculated from the 1970 Jacchia 
atmosphere model described in reference (10) . The following is a list 
of the effects causing atmospheric density variations used in the model; 



(a) variations with an 11.5 year solar cycle 

(b) variations with daily changes in solar activity 

(c) diurnal variations 

(d) variations with geomagnetic activity 

(e) semiannual variations 

(f) seasonal- latitudinal variations of the lower thermosphere 

(g) seasonal- latitudinal variation of helium 

The inputs are the sun and satellite right ascension and declination, 
number of days from January 1, 1970 and vehicle's altitude (km) above 

3 

the surface of the earth. The output is the atmospheric density (kg/m ) 
at the altitude of the satellite. The model calculates atmospheric 
densities for altitudes of 125 km to 700 km with a maximum error of 5% 
when compared to tabulated density values. 
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N=8 K-N -K , K . ^ , K x,. , K , . . 

Pn 

N»1 K=0 e 

I 


N=8 K=N R N+2 K K K 

2 Z -(N+l)(|-)”^^ [g^ cosKX + h^ sinKX] p^ (cos6 ) 
N=1 K=0 e IN n w « n 


Bq, e,, and e are the magnetic spherical force components in a 
0 A r 

geocentric coordinate reference with respect to the Greenwich time line 

as shovm in Figure 9. X„ is the east longitude from the Greenwich line 

M 

and 0.. is the co latitude. R is the satellite orbit radius and R the 
M e 

the earth radius. 

The function p^^ (cos (6^)) is defined as 


1 ^ e^(N-K)!(l-v-^)-^ 'l/2 ,(K4-N) 

‘’n^ ’ (N+K): ‘ 


where 


V = cos 6, 


e„ = 1 for K = 0 
= 2 for K > 1 

K. 


The epoch time for this model is 1965. 


Since the magnetic 


K K 

field varies with time, the coefficients and h^^ are functions of time 


^ (t) = (t^) + (- 

K K 

Reference 12 lists the magnetic coefficients g^ and h^j at epoch 1965 

•K 

with associated secular coefficients g^ and h^. 
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The earth's magnetic field components can be transformed to a 
geocentric Greenwich coordinate system with Z passing through the 

(j 

Greenwich line as shown in Figure 9 by 
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G->S 


cos X, 
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sin6 


-sin X 
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(50) 


-cos6„ sinX^, 
M M 


sin6.. sinX^, 
M M 
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-COS0.. cosX„ 
M M 


COS0 
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sin0^, cosX„ 
M M 


(51) 


Tlie Greenwich geocentric coordinate system is related to the 
geocentric inertial coordinate system by 
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(52) 


where 


cos X 


-sinX 


[L] 


i-k; 


(53) 


sin X 


cosX 
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X is defined as A, + t. A. is the initial angle of the 
o 1 e X 

Greenwich line from the vernal equinox, is the earth rotation rate. 
The magnetic torque is computed from 

M = M X B (54) 

m 

where M and B are the satellite magnetic moment and the earth magnetic 
field vector, respectively. 

Solar Radiation Pressure 

Solar radiation torques are due to the incoming solar radiation 
flux of the sun. The solar torques are functions of the distance from 
the sun, satellite surface geometry, and surface reflectivity. Only 
direct radiation from the sun is considered. Earth and atmospheric 
reflected radiation along with the satellite radiation are ignored. In 
this paper the secular and periodic terms are separated by averaging one 
orbit, assuming a constant inertial solar radiation force. 

The physical model used is similar to the model in reference 
(13) . Figure 10 shows the geometry of the model. The solar radiation 
force is due to the reflected and absorbed components of the Incoming 
sun radiation flux 

dF = dF. + dF (55-A) 

s X r 

where F. and F are the absorbed and reflected solar radiation forces, 

X r 

respectively. Assuming secular reflection, angles ^hd ere equal. 

Yi = Y 2 = Y (55-B) 
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From Figure 10, Is in the opposite direction of the incoming radiation 
flux vector, ii is defined as the unit vector normal to the surface 
element dA. The Incident force component can be written as 


dfj - - (S^ + S^) Pg dA 


(56) 


where 

S. = satellite surface absorption coefficient 

A 

S = satellite surface reflection coefficient 
r 

P = solar radiation pressure 
dA = element area 

The reflected force component is 


dF = - S P^dA e. 
r r S 2 


(57) 


where 62 is defined in Figure 10. From equations 56 and 57 the total 
solar radiation force is 


dF^ = - ^ V "1 + \ 


(58) 


The total solar radiation force can be expressed in terms of the 
incoming flux vector (e^) and the surface normal vector (n) . 

dF = - P_dA [(S + 2S ) COSY n+S.n + (e x ii) ] (59) 

S Q IT xV X 

The solar force components expressed in geocentric inertial 
components are 
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( 60 ) 
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Since the Initial secular trends were of interest, the periodic and 
secular terms were separated whenever possible. For the solar 
radiation force the inertial force components were assumed constant 
over one prbit. The inertial components x^ere transformed to the 
orbital reference coordinate system and averaged over an orbit as a 
function of true anomaly. 

F 

avg 2 tt 

0 


F da 
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( 61 ) 
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where 




avg 


(cosa cosi cosg - sina sing) 
sini cosg 

(-sina cosi cosg - cosa sin3) 


0 

-cosi 

0 


(-cosa cosi sing - sina cosg) 
-sini sing 

(sina cosi sing - cosa cosg) 


( 62 ) 


( 63 ) 


( 64 ) 
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Equation (64) assumes i, F , F , and F are constant over one orbit. 

X y z 

The ascending node was not assumed constant. As the orbit inclination 
approaches 90 degrees, 6 approaches zero. Therefore, if 0 is assumed 
constant the transformation, equation (64), after averaging becomes 


'Vii 


avg 


0 

sini COS0 
0 


0 0 

-cosi -sini sin0 

0 0 


(65) 


To predict accurately long term effects the apparent motion of the 
sum needs to be included. In this study we were concerned with initial 
secular effects in determining the trends of the satellite attitude 
motion. The apparent motion of the sum was assumed constant. 

After averaging the solar forces, the orbital force components 
were transformed to the body reference coordinates. The final averaged 
force components, therefore, can be expressed in matrix notation as 


r- 1 



F i 
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F 

X 

F 
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avg 

F 

z 


b 

L 


The solar torques are computed from the cross-product 

N 

dM_ =2 (p, X dF ) (67) 

X=1 

where is the distance between the center of mass and center of 
pressure of the satellite. N is the number of satellite components. 

The shadow model is used to determine the points at which the 
satellite enters and leaves the earth shadow by assuming the projection 




$5 

of the shadow to be a cylinder. The model is shown in Figure 11 and 
asstnnes the earth, sun, and satellite to be coplanar. 

L_ is the unit vector from the earth to the sun. 

O 

L_ = sin0 cosi i_ + sin0 sini + cos0 (68) 

S s si s s'^I si 

From Figure 11 the angle D can be defined as 

cos D - — (69) 

For cos D > 0, the satellite is in sunlight; for cos D < 0 the satellite 
is in the shadow of the earth. Also from Figure 11 angle E can be 
written as 

sin E = R /r (70) 

e 

Considering Figure (11) geometry the following condition exists; 
if (D + E) < 180°; satellite is in sunlight 
if (D + E) > 180°; satellite is in earth's shadow 

The identity 

sin (D + E) * sin D cos E + sin E cos D ' (71) 

was used in the computer program for this test. 








CHAPTER V 

APPLICATION TO SKYLAB 



The Euler moment equations discussed in Chapter III were solved 
numerically using an IBM 370/168 computer for the Skylab spacecraft. 

The equations were solved using a fourth-order Runge-Kutta method and 

a fourth-order modified Adaras-Bashforth predictor corrector method with 

(14) ^ 

a constant integration step size. The fourth order Runge-Kutta 

method was used to compute the starting values for the Adams-Bashf orth 

predictor-corrector. After the starting values were computed the 

predictor-corrector algorithim was used to integrate the equations of 

motion. 

The Euler's moment equations were solved for the slow tumbling 
Skylab spacecraft using initial conditions in Appendix A. The slow 
tumble case simulates a slow spin approximately about the major axis. 

The slow tumble was solved for the torque-free motion, gravity-gradient, 
aerodynamic drag, solar radiation pressure, and gravity-gradient with 
aerody o torque. Results from each perturbation are discussed in 
Chaptei VI. ^ 

Since the solution of Euler's attitude equations used excessive 
computer time only the initial effect of the perturbations on the 
spacecraft were simulated. The torque-free motion is discussed first, 

"x 

followed by the environmental effects. 

The magnetic effects were not studied since data was not 

1 

' % 

available on the magnetic residual moments of Skylab; they were assumed 

to be negligible. Solar radiation effects were simulated and found to I 

be negligible compared to aerodynamic drag and gravity -gradient torques. 




CHAPTER VI 
SIMULATION RESULTS 



Torque-Free Results 

Since the general torque-free solution to Euler's attitude 
equations are elliptic functions, equation (12) was solved numerically 
to determine the force-free motion in the orbital and satellite 
coordinate system for Skylab, using the initial conditions in Appendix 
A. The motion in the orbital coordinate system of the satellite 
angular momentum vector is represented in Figure 12 for approximately 
five orbits. The nutation and precession angles are defined by 
equation (18). The rapid motion at 40, 90, and 140 minutes is caused 

by the singularity of the Euler Angle 9. As 0 approaches 90 degrees 

• • 

the Euler rates (<|),i/j) become larger as shovm by equation (13). This 
causes the angular momentum to be transferred from either x or y to the 
y or X axis inducing a large angular momentum precession rate for a 
short duration. At t = 180 minutes the precession angle changes 
direction and the angular momentum vector passes above the orbit plane 
(nutation angle > 90°), The nutaltibn angle varies from 70 to 116 
degrees and the precession angle varies from 40 to 326 degrees. Figure 
12 illustrates the motion of the angular momentum vector in the orbital 
coordinate system. 

The nadir angle (A) represents the angular mr^tion of the minor (x) 
axis with respect to the orbit radius vector. Figure 13 illustrates the 
nadir angle for the slow tumble torque free solution. From Figure 13 
the nadir angle oscillates between 0 and 180 degrees. Tlie nadir angle 





Nadir Angles (degrees) 





Figure 13/ Nadir Angle for Torque-Free Solution 
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period is approximately two orbits. The angular momentum and energy 

2 2 2 
for the slow tumble case were 4227 kg»m /s and 2.37 kg*m /s • 

Gravity-Gradient Results 

The result of the gravity-gradient torque on the nadir angle is 
illustrated in Figure 14 which shows that the nadir angle is bound 
between 0 and 40 degrees. Comparing Figures 13 and 14 the gravity- 
gradient torque tends to align the minor axis along the orbit radius 
vector. From inspection of equation (33) an equilibrium position (zero 
torque) state exists when the minor axis is aligned along the orbit 
radius vector. This state occurs when the nadir angle is 0 or 180 
degrees (6 = ±90°). Another equilibrium position exists when all 
principal axes are aligned along the orbital coordinate system axes. 

As shown by equation (33) any misalignment of the principal axes from 
the orbital coordinate system induces a torque. Once a misalignment 
occurs gravity- gradient torques will attempt to orient the satellite 
toward an equilibrium state. Comparing the torque-free solution in 
Figure 13 and the gravity-gradient solution in Figure 14, the amplitude 
of the nadir angle oscillation is reduced from 180 degrees for the 
torque- free case to 40 degrees for the gravity-gradient case. The nadir 
angle period of oscillation decreases from 184 minutes to 66 minutes. 
Thus, by bounding the amplitude of the nadir angle oscillations gravity- 
gradient torque causes the frequency of the oscillations to increase. 

Figure 15 illustrates the gravity-gradient effect on Skylab’s 
energy and angular momentum. Gravity-gradients cause the energy and 
angular momentum to become cyclic with a period of 34.1 minutes or 
approximately .36 of an orbital period. From the initial conditions the 
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spacecraft was near a zero-torque (9 = -80°) position initially. As 
the spacecraft moves from the initial position gravity-gradient torques 
attempt to restore it, decreasing the angular momentum and energy. As 
a result angular momentum and energy are maximum when the satellite 
nadir angle is minimum and vice versa. Figures 14 and 15 illustrate 
this effect. 

The motion of the angular momentum vector in the orbital coordinate 
system is shown in Figure 16. The gravity-gradient torque causes the 
precession angle and rate to become periodic with an approximate period 
of fifty minutes. The singularity due to Euler angle 0=0 degrees in 
the torque-free case no longer occurs since the Euler angle 0 is bound 
bet^^een -45 and -90 degrees. The precession angle represents the change 
of angular momentum in the orbital coordinate system lii the orbit 
normal- velocity vector plane. Examination of the gravity-gradient 
torque equations in the orbital coordinate system, equations (34), for a 
symmetric satellite about the minor axis (B=C) indicates why the regular 
precession occurs. Making the assumption (B=C) for Sky lab the gravity- 
gradient torque equations in the orbital coordinate system become 

(M^)^ = sin20 sinij; (A-C) 3y/2R^ 

(M ) = sin20 cosi (C-A) 3y/2R^ (72) 

y o 


(M ) = 0 

z o 

and M^ are the torque components in the orbit normal-velocity vector 
plane and are functions of 0 and t|;. Since 9 is bound between -40 and 
-90 degrees the Euler angle i|^ controls the direction of the precession. 

The nutation angle represents the motion of the angular momentum 
vector with respect to the orbit radius vector. Since the z gravity- 
gradient torque component in the orbital coordinate system is zero the 
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the nutation angle reflects the change in the total angular momentum 
magnitude. The nutation angle is shown in Figure 16. 

Gravity-gradient causes the angular momentum and energy to 
become cyclic, attempts to align the minor axis along the orbit radius 
vector, and causes a regular precession of the angular momentum about 
the orbit radius vector. 

Aerodynamic Drag Results 

The aerodynamic torque equations are complex and difficult to 
examine analytically. However, by examining the trends of various 
spacecraft attitude parameters, valuable insight into the problem can 
be gained. 

Figure 17 illustrates the motion of the minor axis with respect to 
the orbit radius vector. As compared with the torque-free case, Figure 
13, aerodynamic torques tend to damp out the oscillations and drive the 
minor axis perpendicular to the orbit radius vector. Since the free 
streamline velocity vector is approximately parallel to the orbit velo- 
city vector, aerodynamic drag would orient the spacecraft toward a 
position of minimum resistance. From Figure 17 this position would 
appear to occur when the minor axis is in the orbit normal velocity 
vector plane. This can be illustrated by examining the nadir angle. 

The nadir angle amplitude continues to decrease in magnitude until 310 
minutes and is oscillating about the velocity vector tangent (nadir angle 
of 90°). 

Figure 18 illustrates the change in angular momentum and energy. 

As shown the angular momentum and energy decrease until approximately 
310 minutes. After this time the angular momentum and energy begin to 
Increase. 
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An explanation for the above trends can be found by examining 
the change in angular momentum in the spacecraft body coordinate 
system. Figure 19 illustrates the angular momentum change for the 
major axis (B) and the intermediate axis (C) . The angular momentum 
change for the minor axis (A) was very small. As shovm in Figure 19 
initially the spacecraft was principally spinning about the major axis. 
As time progressed the aerodynamic effect caused the angular momentum 
to shift from the major axis to the intermediate axis until the space- 
craft was essentially spinning about the intermediate axis at approxi- 
mately 310 minutes. From attitude dynamics a spacecraft spinning 
about the intermediate axis with a perturbation is in an unstable 
mode. Reference 15 illustrates why an unstable mode exists for motion 
about the intermediate axis with a perturbation. After reaching this 
unstable mode the spacecraft quickly changes its momentum state and 
within an orbit is spinning about the major axis. It is interesting to 
note that the spacecraft is spinning about the major axis in the 
direction opposite to its initial one. This change .in angular momentum 
state causes the nadir angle, angular momentum, and energy to increase. 

The aerodynamic drag effect on the position of the angular 
momentum vector in the orbital coordinate system is illustrated in 
Figure 20. As compared to the torque-free solution, the oscillatory 
motion of the nutation angle is damped out and the momentum vector tends 
to align along the radius vector. The nutation angle continues to 
decrease until approximately 310 minutes. As the minor axis is being 
forced to become parallel to the orbit’s velocity vector tangent the 
momentum is being transferred from the major to intermediate axis as 
shown in Figure 19 and explains the reason for the decrease. The 


Time (minutes) 


Figure 19. Effect of Aerodynamic Drag on the Y and Z Angular 
Momentum Compqnents in the Body Coordinate System. 






Time (minutes) 

Figure 20. Aerodynamic Drag Effect on the Precession and Nutation Angles 
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interaediate axis is oriented in the general direction of the orbit's 
radius vector. As a result the nutation angle decreases until the 
angular momentum is transferred back to the principal axis. As the 
angular momentum is transferred back to the major axis the nutation 
angle increases. The aerodynamic effect on the precession angle is to 
cause it to precess slowly about the orbit radius vector as shown in 
Figure 20. 

Combined Effects Oj. Gravity-Gradient and Aerodynamic Torques 

The combined effect of aerodynamic and gravity-gradient torques on 

the energy and angular momentum is illustrated in Figure 21. Since 

gravity- gradient is larger than aerodynamic drag the cyclic motion due 

to gravity-gradient is predominant. The periodic type motion of the 

angular momentum and energy resulting from gravity-gradient becomes a 

random type oscillation due to the damping effect of aerodynamic drag. 

The aerodynamic torque reduces the amplitude of the gravity- gradient 

induced oscillations. Initially the amplitude of oscillations varied 

2 

between 600 and 4200 kg*m /s; as time progressed the amplitudes varied 

2 

between 2700 and 4200 kg*m /s. The lower limit of the gravity-gradient 
oscillation was affected more severely. This was due to the aerodynamic 
effect of aligning the minor axis along the velocity vector. 

Comparing the nadir angle of gravity-gradient (Figure 14) and 
gravity-gradient with aerodynamic drag (Figure 21) illustrates the 
aerodynamic effect. Aerodynamic drag reduces the amplitude of oscilla- 
tion of the angular momentum and nadir angle which were induced by the 
gravity- gradient torque. This torque tends to align the minor axis 
along the radius vector. The aerodynamic torque attempts to align the 



minor axis along the orbit velocity vector tangent perpendicular to 
the orbit radius. The combined effect of gravity-gradient and 
aerodynamic torque is to cause the minor axis to oscillate between 
the radius vector and orbit velocity vector tangent as illustrated in 
Figure 22. Since the gravity- gradient torque is larger than aerodynamic, 
its effect is predominant. This causes the minor axis to oscillate 
nearer the orbit radius vector between 25 and 30 degrees. 

Since the gravity-gradient torque component in the orbit radius 
vector is zero the change in nutation angle illustrated in Figure 23 is 
due to the aerodynamic torque. The effect of the aerodynamic torque is 
to maintain the angular momentum vector below the orbit velocity vector- 
orbit normal plane. The precession angle as shown in Figure 23 has a 
period of approximately 46 minutes. This is similar to gravity-gradient 
effect shown in Figure 14. The difference is aerodynamic torque 
eliminates the change in the precession angle direction. As a result 
the combined effect of gravity-gradient torque with aerodynamic drag is 
to cause the angular momentum vector to precess regularly about orbit 
radius vector pointing down the radius vector. 


Nadir Angles (degrees) 
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Figure 22. Nadir Angle for Aerodynamic Drag and Gravity “Grad lent. 
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Figure 23. Aerodynamic Drag and Gravity-Gradient Effect on the Precession 

and Nutation Angles. 
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CHAPTER VII 

CONCLUSIONS AND RECOMMENDATIONS 

An investigation of environmental perturbations for an asymmetric 
slow tumbling satellite has been presented. Environmental perturbation 
sources considered were gravity-gradient, aerodynamic drag, solar 
radiation pressure and magnetic field interactions. For the Skylab 
spacecraft results assuming gravity-gradient and aerodynamic torques 
were presented. Solar radiation pressure and magnetic torques were 
small and neglected. 

Gravity- grad lent torque causes the nadir angle to become bounded 
about the orbit radius vector and causes the energy, angular momentum, 
and precession angle to become cyclic. Aerodynamic drag initially 
decreased the angular momentum, energy, and nadir angle and drives the 
spacecraft into an unstable mode by transferring angular momentum from 
the major axis to the intermediate axis. Aerodynamic drag also induces 
slow precession rate about the orbit radius vector and causes the 
nutation angle to decrease. The combined effects of gravity-gradient 
with aerodynamic torques were to cause the gravity-gradient induced 
oscillation amplitudes to decrease and the periodic motion to become 
random type oscillations. Gravity-gradient torques were greater than 
aerodynamic drag and its effect was predominant. 

For long term predictions of satellite tumbling motion, analytical 
techniques must be investigated. Computer simulation for this type of 
motion becomes very expensive since small time intervals are required to 
solve the nonlinear equations of motion. Future models of tumbling 
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spacecraft should include energy dissipation models to account for 
internal energy dissipation. 
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APPENDIX A 

Initial Conditions for the Slow Tumble Mode 
A. Slow Tumble Orientation 

The Skylab spacecraft was assumed to be in the following attitude 
jjj = 1.07 degrees 
9 = -79.96 degrees 
(j) = 12.85 degrees 

This was the attitude Skylab was believed to be in when it became 

passive. For the slow spin case a negative angular velocity equal in 

magnitude to the orbital rate about the orbit normal was assumed. 

0) = +.00112256 radians /second 

o 

O) = to j 
o o o 

Using the transformations and Figure (4) the required initial 

Euler rates can be computed. 

- 0.0 

0 = 03 cos 
o 

(f) = 03 sin i|3 Cos 0 
o 

Assuming 03^ is the negative of the orbital rate the initial Euler rates 
for the flow tumble case are 

if = 0,0 

9 = -.0011229 radians/second 

* -6 

(}> = -3.656 X 10 radians/second 


V OF 
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APPENDIX B 

Skylab’s Aerodynamic Drag Moment Coefficients 



^0 


^2 

^3 

"4 

% 

% 


•>2 

^3 

Cy (pitching moment 

coefficient) 

-1.289 

0.135 

-0.195 

-0.009 

0.080 

0.007 

-0.037 

-0.004 

0.009 

-0.011 

0.022 

-0.005 

0.011 

0.442 

-0.073 

0.084 

0.000 

-0.076 

-0.011 

0.003 

0.009 

-0.022 

0.015 

-0.036 

0.007 

-0.020 


-0.037 

-5.146 

0.011 

0.350 

0.01^ 

-0.010 

-0.007 

0.199 

-0.017 

0.120 

0.002 

-0.06 

0.009 


-0.044 

0.779 

-0.007 

0.123 

0.005 

-0.032 

-0.001 

-0.049 

-0.012 

-0.021 

0.001 

0.032 

-0.001 






c^ (yawing 

moment 

coefficient) 

-0.958 

-0.067 

-0.266 

0.019 

0.016 

-0.006 

0.012 

-0.096 

-0.055 

-0.026 

-0.054 

0.018 

-0.009 


0.269 

0.006 

0.207 

-0.007 

-0.032 

0.006 

-0.014 

0.055 

0.031 

0.014 

0.068 

-0.027 

0.012 


0.166 

-0.140 

-0.014 

0.280 

-0.042 

-0.074 

0.046 

4.210 

0.022 

0.454 

-0.052 

-0.055 

0.003 


0.155 

0.072 

0.003 

-0.069 

-0.015 

0.032 

0,004 

-0.467 

0.004 

-0.134 

-0.010 

0.003 

-0.001 




(rolling moment coefficient) 


-0.056 

-0.071 

-0.038 

0.018 

-0.035 

0.013 

0.078 

0.047 

0.058 

-0.021 

0.042 

-0.016 

0.053 

1.136 

0.024 

0.943 

-0.016 

0.430 

0.025 

-0.210 

0.009 

-0.162 

0.005 

-0.099 


0.009 

-0,012 

-0.001 

-0,017 

0,002 

-0.011 

-0.009 

0,018 

0.005 

0.028 

-0,000 

0.018 

0.250 

-0.005 

-0.039 

0.005 

-0.052 

-0.005 

-0.055 

-0.005 

0.019 

0.000 

0.007 

-0.002 


0.006 

-0.016 

-0.007 


-0.003 
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APPENDIX C 

Skylab's Orbit and Satellite Parameters 

Orbital Parameters 

Orbit Inclination = 50 degrees 

Eccentricity = 0 

Altitude = 435,5 kilometers 

Right Ascension of Ascending Node = 233.2 degrees 
Orbital Period = 93.23 minutes 

5 3 

Earth's gravitational constant = 3.986 S 10 kilometer / 

second‘d 

Earth's mean radius - 6378.0 kilometers 

Satellite Parameters 

Principal Moment of Inertias 

I = 7.93321 X 10^ kilograms - meter 2 

XX 

I = 3.767828 x 10 kilograms •- meter 2 

yy 

I = 3.694680 X 10 kilograms - m.eter 2 
zz 

Aerodynamic Reference Area 

2 

Surface Reference Area = 79.46 meter 


Reference Diameter = 10.058 meter 


